#setwd('/afs/inf.ed.ac.uk/user/s17/s1725186/Documents/PhD-Models/FirstPUModel/RMarkdowns')
library(tidyverse) ; library(reshape2) ; library(glue) ; library(plotly) ; library(dendextend)
library(RColorBrewer) ; library(viridis) ; require(gridExtra) ; library(GGally) ; library(gplots)
library(WGCNA)
library(expss)
library(polycor)
library(knitr)
SFARI_colour_hue = function(r) {
pal = c('#FF7631','#FFB100','#E8E328','#8CC83F','#62CCA6','#59B9C9','#b3b3b3','#808080','gray','#d9d9d9')[r]
}
clustering_selected = 'DynamicHybridMergedSmall'
#clustering_selected = 'DynamicHybrid'
cat(paste0('Using clustering ', clustering_selected))
## Using clustering DynamicHybridMergedSmall
Load preprocessed dataset (preprocessing code in 19_10_14_data_preprocessing.Rmd) and clustering (pipeline in 19_10_15_WGCNA.Rmd)
# Gandal dataset
load('./../Data/Gandal/preprocessed_data.RData')
datExpr = datExpr %>% data.frame
DE_info = DE_info %>% data.frame
# GO Neuronal annotations
GO_annotations = read.csv('./../../../PhD-InitialExperiments/FirstYearReview/Data/GO_annotations/genes_GO_annotations.csv')
GO_neuronal = GO_annotations %>% filter(grepl('neuron', go_term)) %>%
mutate('ID'=as.character(ensembl_gene_id)) %>%
dplyr::select(-ensembl_gene_id) %>% distinct(ID) %>%
mutate('Neuronal'=1)
# SFARI Genes
SFARI_genes = read_csv('./../Data/SFARI/SFARI_genes_08-29-2019_with_ensembl_IDs.csv')
SFARI_genes = SFARI_genes[!duplicated(SFARI_genes$ID) & !is.na(SFARI_genes$ID),]
# Update DE_info with SFARI and Neuronal information
genes_info = DE_info %>% mutate('ID'=rownames(.)) %>% left_join(SFARI_genes, by='ID') %>%
mutate(`gene-score`=ifelse(is.na(`gene-score`), 'None', `gene-score`)) %>%
left_join(GO_neuronal, by='ID') %>%
mutate(Neuronal=ifelse(is.na(Neuronal), 0, Neuronal)) %>%
mutate(gene.score=ifelse(`gene-score`=='None' & Neuronal==1, 'Neuronal', `gene-score`),
significant=padj<0.05 & !is.na(padj))
# Dataset created with DynamicTreeMerged algorithm
dataset = read.csv(paste0('./../Data/Gandal/dataset_', clustering_selected, '.csv'))
dataset$Module = dataset[,clustering_selected]
# GO DEA
# load('./../Data/Gandal/GO_DE_clusters.RData')
rm(DE_info, GO_annotations, clusterings)
Check if for each gene, the module with the smallest distance is their corresponding module
mm = dataset %>% dplyr::select(starts_with('MM.'), starts_with('MMgray')) %>% dplyr::rename(MM.gray = MMgray)
rownames(mm) = dataset$ID
original_max_membership = gsub('MM.', '', colnames(mm)[max.col(mm, ties.method='first')])
cat(paste0('For ', round(100*mean(original_max_membership == gsub('#','',dataset$Module))),
'% of the genes, their assigned module corresponds to the module with the highest Module Membership'))
## For 53% of the genes, their assigned module corresponds to the module with the highest Module Membership
Apparently this is not the case. Someone asked this in a Bioconductor question and Peter Langfelder answered that this is because WGCNA assigns module labels using dynamic tree cut of hierarchical clustering tree that is based on the Toplogical Overlap Measure. TOM results in similar but not quite the same similarity as correlation, hence for some genes the assigned module may differ from the module with highest kME. But that in all, he doesn’t worry about the module assignment vs. max. kME differences in his own analyses, and he recommends not worrying about it to others as well.
Even when the maximum MM module and the assigned module don’t match, the assigned module is almost always one of the highest ranked Modules by MM, so it’s not that bad that these two metrics don’t match.
mm_rank = matrix(0, nrow=nrow(mm), ncol=ncol(mm))
for(i in 1:nrow(mm_rank)) mm_rank[i,] = rank(-mm[i,], ties.method='min')
colnames(mm_rank) = colnames(mm)
rownames(mm_rank) = rownames(mm)
assigned_module_rank = rep('', nrow(mm_rank))
for(i in 1:nrow(mm_rank)){
assigned_module_rank[i] = mm_rank[i, paste0('MM.',gsub('#','',dataset$Module[i]))]
}
# Remove genes corresponding to the gray module
assigned_module_rank = assigned_module_rank[dataset$Module!='gray']
plot_data = assigned_module_rank %>% table %>% data.frame
colnames(plot_data) = c('rank', 'frequency')
plot_data = plot_data %>% mutate(rank = as.numeric(as.character(rank))) %>%
arrange(rank) %>% mutate('percentage'=round(100*frequency/sum(frequency),2))
ggplotly(plot_data %>% ggplot(aes(rank, frequency, fill=percentage)) + geom_bar(stat='identity') +
theme_minimal() + theme(legend.position='none') + ggtitle('Ranking of assigned module by MM') +
xlab('MM rank by assigned Module'))
The palette on the top of the heatmap represents the size of the module, the darker the colour, the larger the module
module_size = dataset %>% group_by(Module) %>% tally %>% mutate(Module = gsub('#', '', Module))
module_size$quant = cut(module_size$n, breaks=9, labels=FALSE)
module_size$meanExpr = sapply(module_size$Module, function(m){mean(rowMeans(datExpr)[dataset$Module==m | dataset$Module==paste0('#',m)])})
module_size$quantME = cut(module_size$meanExpr, breaks=9, labels=FALSE)
heatmap.2(table(original_max_membership, gsub('#','',dataset$Module)), symm=TRUE, dendrogram='none', keysize=1,
trace='none', scale='none', col=brewer.pal(9,'YlGnBu'), xlab='Assigned Module', ylab='Highest MM',
ColSideColors = brewer.pal(9,'PuBu')[module_size$quant])
Scaling by rows it’s easier to see that genes assigned to the larger modules sometimes have a higher module membership with smaller modules
heatmap.2(table(original_max_membership, gsub('#','',dataset$Module)), dendrogram='none', keysize=1, symm=TRUE,
trace='none', scale='row', col=brewer.pal(9,'YlGnBu'), xlab='Assigned Module', ylab='Highest MM',
ColSideColors=brewer.pal(9,'PuBu')[module_size$quant])
Module Membership by gene (selecting a random sample so it’s not that heavy)
The membership of some module seems to be related to the level of expression of the genes
set.seed(123)
plot_mm = mm %>% sample_frac(0.05) %>% as.matrix
colnames(plot_mm) = gsub('MM.','', colnames(plot_mm))
plot_mm = plot_mm[order(rowMeans(plot_mm)), module_size$Module[order(module_size$n)]]
heatmap.2(plot_mm, xlab('Modules ordered by Size'), ylab('Genes ordered by mean expression'), keysize=1,
ColSideColors = brewer.pal(9,'PuBu')[module_size$quant[order(module_size$n)]],
RowSideColors = brewer.pal(9, 'RdPu')[cut(rowMeans(plot_mm), breaks=9, labels=FALSE)],
trace='none', dendrogram='none', scale='none', Rowv=FALSE, Colv=FALSE, col=brewer.pal(9,'Spectral'))
Letting the heatmap order the genes and modules by distance, it seems like the module size is not an important factor but the mean expression is. The modules in the right leg of the dendrogram seem to be associated to genes with low expression and the ones on the left leg to genes with high expression
heatmap.2(plot_mm, trace='none', col=brewer.pal(9,'Spectral'), dendrogram='column',
RowSideColors = brewer.pal(9, 'RdPu')[cut(rowMeans(plot_mm), breaks=9, labels=FALSE)],
ColSideColors = brewer.pal(9,'PuBu')[module_size$quant[order(module_size$n)]])
The biggest modules have the most extreme memberships (both positive and negative)
plot_mm = plot_mm[order(rowMeans(abs(plot_mm))), module_size$Module[order(module_size$n)]]
heatmap.2(abs(plot_mm), ColSideColors = brewer.pal(9,'PuBu')[module_size$quant[order(module_size$n)]], keysize=1,
RowSideColors = brewer.pal(9, 'RdPu')[cut(rowMeans(plot_mm), breaks=9, labels=FALSE)],
trace='none', dendrogram='none', scale='none', Rowv=FALSE, Colv=FALSE, col=brewer.pal(9,'YlOrRd'))
Checking if the mean expression of the modules plays a role
plot_mm = plot_mm[order(rowMeans(plot_mm)), module_size$Module[order(module_size$meanExpr)]]
heatmap.2(plot_mm, ColSideColors = brewer.pal(9,'RdPu')[module_size$quantME[order(module_size$meanExpr)]], keysize=1,
RowSideColors = brewer.pal(9, 'RdPu')[cut(rowMeans(plot_mm), breaks=9, labels=FALSE)],
trace='none', dendrogram='none', scale='none', Rowv=FALSE, Colv=FALSE, col=brewer.pal(9,'Spectral'))
plot_mm = plot_mm[order(rowMeans(plot_mm)), module_size$Module[order(module_size$meanExpr)]]
heatmap.2(plot_mm, ColSideColors = brewer.pal(9,'RdPu')[module_size$quantME[order(module_size$meanExpr)]], keysize=1,
RowSideColors = brewer.pal(9, 'RdPu')[cut(rowMeans(plot_mm), breaks=9, labels=FALSE)],
trace='none', dendrogram='none', scale='none', Rowv=TRUE, Colv=FALSE, col=brewer.pal(9,'Spectral'))
The mean expression pattern seems to be clearer in the samples than in the modules
The colors in the histogram seem to have inverted, with the red being the highest MM and blue de lowest (I have no idea how this happened)
plot_mm = plot_mm[order(rowMeans(plot_mm)), module_size$Module[order(module_size$meanExpr)]]
heatmap.2(plot_mm, ColSideColors = brewer.pal(9,'RdPu')[module_size$quantME[order(module_size$meanExpr)]], keysize=1,
RowSideColors = brewer.pal(9, 'RdPu')[cut(rowMeans(plot_mm), breaks=9, labels=FALSE)],
trace='none', dendrogram='column', scale='none', col=brewer.pal(9,'Spectral'))
Because of the weird inversion in the heatmap palette from above, check if there’s a relation between the gene’s level of expression and the average level of expression of the module it is assigned to: there is and it is positive, so something weird is happening in the heatmaps above
I’m not sure how to compare these two plots …
plot_data = data.frame('ID' = dataset$ID, 'GeneMeanExpr'=rowMeans(datExpr), 'Module' = gsub('#','',dataset$Module)) %>%
left_join(module_size, by='Module') %>% mutate('ModuleMeanExpr' = meanExpr)
MM_module_size = data.frame('Module' = original_max_membership) %>% group_by(Module) %>% tally()
MM_module_size$meanExpr = sapply(MM_module_size$Module, function(m){mean(rowMeans(datExpr)[original_max_membership==m |
original_max_membership==paste0('#',m)])})
MM_plot_data = data.frame('ID' = dataset$ID, 'GeneMeanExpr'=rowMeans(datExpr), 'Module' = gsub('#','',dataset$Module)) %>%
left_join(MM_module_size, by='Module') %>% mutate('ModuleMeanExpr' = meanExpr)
p1 = plot_data %>% ggplot(aes(ModuleMeanExpr, GeneMeanExpr)) + geom_point(alpha=0.2, color=gsub('#gray','gray',paste0('#',plot_data$Module))) +
geom_smooth(color='gray') + theme_minimal() + xlab('Assigned Module Mean Expression') + ylab('Gene Mean Expression') +
ggtitle(paste0('Cor = ', round(cor(plot_data$ModuleMeanExpr, plot_data$GeneMeanExpr),3),
' Regression slope = ', round(lm(GeneMeanExpr ~ ModuleMeanExpr, data=plot_data)$coefficients[[2]],2),
' R^2 = ', round(cor(plot_data$ModuleMeanExpr, plot_data$GeneMeanExpr)^2,3)))
p2 = MM_plot_data %>% ggplot(aes(ModuleMeanExpr, GeneMeanExpr)) + geom_point(alpha=0.2, color=gsub('#gray','gray',paste0('#',plot_data$Module))) +
geom_smooth(color='gray') + theme_minimal() + xlab('Highest MM Module Mean Expression') + ylab('Gene Mean Expression') +
ggtitle(paste0('Cor = ', round(cor(MM_plot_data$ModuleMeanExpr, MM_plot_data$GeneMeanExpr),3),
' Regression slope = ', round(lm(GeneMeanExpr ~ ModuleMeanExpr, data=MM_plot_data)$coefficients[[2]],2),
' R^2 = ', round(cor(MM_plot_data$ModuleMeanExpr, MM_plot_data$GeneMeanExpr)^2,3)))
grid.arrange(p1, p2, nrow = 1)
rm(p1, p2)
Assigning genes by MM increases the mean expression of the largest modules and decreases the mean expression of the smallest ones
plot_data = module_size %>% dplyr::select(Module, meanExpr, n) %>% rename(AssignedModuleME=meanExpr) %>%
left_join(MM_module_size, by='Module') %>% rename(MMModuleME=meanExpr)
ggplotly(plot_data %>% ggplot(aes(AssignedModuleME, MMModuleME, size=n.x)) + geom_abline(slope=1, intercept=0, color='gray') +
geom_point(color=gsub('#gray','gray',paste0('#',plot_data$Module)), alpha=0.5, aes(id=Module)) + ggtitle('Mean Expression') +
xlab('Mean Expression of Assigned Module') + ylab('Mean Expression of Module with Higest MM') + theme_minimal() + coord_fixed())
There doesn’t seem to be a relation between module size and level of expression, so the result from above probably has no effect on the level of expression patterns found before
plot_data = module_size %>% dplyr::select(Module, meanExpr, n)
ggplotly(plot_data %>% ggplot(aes(meanExpr, n)) + geom_point(color=gsub('#gray','gray',paste0('#',plot_data$Module))) +
geom_smooth(color='gray', alpha=0.1) + theme_minimal())
To see if the strength of the relation between level of expression and assigned module (using the regular WGCNA assignment vs using the maximum MM), I’m going to fit a linear regression to these two module classifications and see which has a better performance.
It seems like the maxMM Module assignment has a stronger relation with mean expression than the originally asigned Modules … but I’m not sure how reliable this results are
lm_data = data.frame('MeanExpr' = log2(rowMeans(datExpr)+1), 'MMModule' = as.factor(original_max_membership),
'assignedModule' = as.factor(gsub('#','',dataset$Module)))
lm_MM = lm(MeanExpr ~ MMModule, data=lm_data[lm_data$MMModule!='gray',])
lm_aM = lm(MeanExpr ~ assignedModule, data=lm_data[lm_data$assignedModule!='gray',])
cat('Results using the maxMM Module:')
## Results using the maxMM Module:
summary(lm_MM)
##
## Call:
## lm(formula = MeanExpr ~ MMModule, data = lm_data[lm_data$MMModule !=
## "gray", ])
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.83981 -0.20169 0.03769 0.21803 0.99391
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.090479 0.014779 209.119 < 2e-16 ***
## MMModule00B0F6 0.011800 0.023610 0.500 0.617244
## MMModule00B6EB 0.057810 0.023325 2.478 0.013205 *
## MMModule00B70C 0.213295 0.022722 9.387 < 2e-16 ***
## MMModule00BB45 0.058565 0.020231 2.895 0.003798 **
## MMModule00BBDD 0.140225 0.022637 6.194 5.99e-10 ***
## MMModule00BD64 0.032734 0.018983 1.724 0.084663 .
## MMModule00BECD -0.074431 0.021525 -3.458 0.000546 ***
## MMModule00BF7D -0.116130 0.021829 -5.320 1.05e-07 ***
## MMModule00C094 -0.186260 0.019680 -9.465 < 2e-16 ***
## MMModule00C0BB 0.271495 0.019518 13.910 < 2e-16 ***
## MMModule00C1A8 0.092988 0.025936 3.585 0.000338 ***
## MMModule4B9FFF 0.057754 0.019862 2.908 0.003645 **
## MMModule53B400 0.179547 0.017727 10.128 < 2e-16 ***
## MMModule75AF00 0.101815 0.022575 4.510 6.52e-06 ***
## MMModule8195FF 0.172702 0.018685 9.243 < 2e-16 ***
## MMModule8EAB00 -0.055586 0.020964 -2.652 0.008020 **
## MMModuleA3A500 -0.080720 0.022258 -3.627 0.000288 ***
## MMModuleA58AFF -0.030488 0.021068 -1.447 0.147887
## MMModuleB4A000 -0.064662 0.024838 -2.603 0.009239 **
## MMModuleC17FFF 0.221578 0.017935 12.355 < 2e-16 ***
## MMModuleC49A00 0.108000 0.018660 5.788 7.27e-09 ***
## MMModuleD29300 0.035031 0.018920 1.852 0.064115 .
## MMModuleD774FD -0.131048 0.024665 -5.313 1.09e-07 ***
## MMModuleDE8C00 -0.077895 0.027104 -2.874 0.004060 **
## MMModuleE76BF3 0.027031 0.018741 1.442 0.149233
## MMModuleE88523 -0.105360 0.023887 -4.411 1.04e-05 ***
## MMModuleF17E4F 0.076575 0.018966 4.037 5.43e-05 ***
## MMModuleF365E6 0.135140 0.018144 7.448 9.95e-14 ***
## MMModuleF8766D 0.015439 0.022055 0.700 0.483925
## MMModuleFB61D7 0.128640 0.021555 5.968 2.45e-09 ***
## MMModuleFD6F86 0.114806 0.019423 5.911 3.47e-09 ***
## MMModuleFF61C5 -0.133784 0.020241 -6.610 3.97e-11 ***
## MMModuleFF64B2 -0.008628 0.019888 -0.434 0.664395
## MMModuleFF699D 0.163659 0.022966 7.126 1.08e-12 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3014 on 16126 degrees of freedom
## Multiple R-squared: 0.1232, Adjusted R-squared: 0.1214
## F-statistic: 66.67 on 34 and 16126 DF, p-value: < 2.2e-16
cat('Results using the assigned Module:')
## Results using the assigned Module:
summary(lm_aM)
##
## Call:
## lm(formula = MeanExpr ~ assignedModule, data = lm_data[lm_data$assignedModule !=
## "gray", ])
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.88606 -0.20640 0.03891 0.22270 0.94648
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.110595 0.026810 116.025 < 2e-16 ***
## assignedModule00B0F6 0.149851 0.041351 3.624 0.000291 ***
## assignedModule00B6EB 0.110688 0.041478 2.669 0.007625 **
## assignedModule00B70C 0.190191 0.030889 6.157 7.58e-10 ***
## assignedModule00BB45 0.043149 0.030760 1.403 0.160713
## assignedModule00BBDD 0.208979 0.043209 4.836 1.33e-06 ***
## assignedModule00BD64 -0.076429 0.027916 -2.738 0.006192 **
## assignedModule00BECD -0.081527 0.041103 -1.983 0.047332 *
## assignedModule00BF7D -0.093901 0.031729 -2.959 0.003086 **
## assignedModule00C094 -0.191956 0.029751 -6.452 1.13e-10 ***
## assignedModule00C0BB 0.221029 0.029092 7.598 3.18e-14 ***
## assignedModule00C1A8 0.331791 0.060508 5.483 4.23e-08 ***
## assignedModule4B9FFF 0.007159 0.029707 0.241 0.809581
## assignedModule53B400 0.085401 0.027705 3.083 0.002056 **
## assignedModule75AF00 0.146637 0.038367 3.822 0.000133 ***
## assignedModule8195FF 0.064730 0.028788 2.249 0.024558 *
## assignedModule8EAB00 -0.109759 0.038774 -2.831 0.004650 **
## assignedModuleA3A500 -0.103021 0.040864 -2.521 0.011709 *
## assignedModuleA58AFF 0.057693 0.039303 1.468 0.142150
## assignedModuleB4A000 -0.039886 0.053020 -0.752 0.451900
## assignedModuleC17FFF 0.116542 0.028091 4.149 3.36e-05 ***
## assignedModuleC49A00 0.051890 0.028551 1.817 0.069162 .
## assignedModuleD29300 0.083385 0.031916 2.613 0.008993 **
## assignedModuleD774FD -0.088074 0.041103 -2.143 0.032147 *
## assignedModuleDE8C00 0.063322 0.053020 1.194 0.232378
## assignedModuleE76BF3 -0.050982 0.028297 -1.802 0.071614 .
## assignedModuleE88523 -0.082744 0.035408 -2.337 0.019459 *
## assignedModuleF17E4F 0.026795 0.028176 0.951 0.341624
## assignedModuleF365E6 0.053883 0.028381 1.899 0.057640 .
## assignedModuleF8766D 0.103814 0.036156 2.871 0.004093 **
## assignedModuleFB61D7 0.252664 0.054412 4.644 3.45e-06 ***
## assignedModuleFD6F86 0.094489 0.028613 3.302 0.000961 ***
## assignedModuleFF61C5 -0.135269 0.033166 -4.079 4.55e-05 ***
## assignedModuleFF64B2 -0.041716 0.038136 -1.094 0.274027
## assignedModuleFF699D 0.176538 0.045017 3.922 8.83e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3069 on 16387 degrees of freedom
## Multiple R-squared: 0.08766, Adjusted R-squared: 0.08577
## F-statistic: 46.31 on 34 and 16387 DF, p-value: < 2.2e-16
The models are quite similar
coefs_MM = summary(lm_MM)$coefficients %>% data.frame %>% mutate('Module' = gsub('MMModule','',rownames(.)))
coefs_aM = summary(lm_aM)$coefficients %>% data.frame %>% mutate('Module' = gsub('assignedModule','',gsub('#','',rownames(.))))
plot_data = coefs_MM %>% left_join(coefs_aM, by='Module') %>% mutate('color'=ifelse(Module=='(Intercept)', 'gray', paste0('#',Module)))
ggplotly(plot_data %>% ggplot(aes(Estimate.x, Estimate.y)) + geom_point(color=plot_data$color, aes(id=Module)) +
geom_abline(slope=1, intercept=0, color='gray') + xlab('Estimate MM Model') + ylab('Estimate assignedModule Model') +
theme_minimal() + ggtitle('Relation between coefficients'))
Assigning genes by MM balances more the size of the modules, increasing the smallest ones and reducing the largest ones. Mean expression doesn’t seem to be a important factor here.
plot_data = module_size %>% dplyr::select(Module, meanExpr, n) %>% rename(AssignedModuleN=n) %>%
left_join(MM_module_size, by='Module') %>% rename(MMModuleN=n)
ggplotly(plot_data %>% ggplot(aes(AssignedModuleN, MMModuleN, size=meanExpr.x)) + geom_abline(slope=1, intercept=0, color='gray') +
geom_smooth(color='gray', se=FALSE) + geom_point(color=gsub('#gray','gray',paste0('#',plot_data$Module)), alpha=0.5, aes(id=Module)) +
ggtitle('Module Size') + xlab('Size of Assigned Module') + ylab('Size of Module with Higest MM') + theme_minimal())
The module assigned to each gene and the module with the highest Module Membership don’t need to be the same and it’s fine, it’s because they are done with different methods
The biggest modules have the strongest MM values (both positive and negative)
There seems to be some relation between the mean level of expression of a gene and MM (this shouldn’t happen…)
If we use maximum MM as a criteria for assigning modules, the resulting modules have more balanced sizes but it looks like the bias by level of expression gets stronger (altgough I’m not sure about how reliable these results are)
sessionInfo()
## R version 3.6.0 (2019-04-26)
## Platform: x86_64-redhat-linux-gnu (64-bit)
## Running under: Scientific Linux 7.6 (Nitrogen)
##
## Matrix products: default
## BLAS/LAPACK: /usr/lib64/R/lib/libRblas.so
##
## locale:
## [1] LC_CTYPE=en_GB.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_GB.UTF-8 LC_COLLATE=en_GB.UTF-8
## [5] LC_MONETARY=en_GB.UTF-8 LC_MESSAGES=en_GB.UTF-8
## [7] LC_PAPER=en_GB.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=en_GB.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] knitr_1.24 polycor_0.7-10 expss_0.10.1
## [4] WGCNA_1.68 fastcluster_1.1.25 dynamicTreeCut_1.63-1
## [7] gplots_3.0.1.1 GGally_1.4.0 gridExtra_2.3
## [10] viridis_0.5.1 viridisLite_0.3.0 RColorBrewer_1.1-2
## [13] dendextend_1.13.2 plotly_4.9.1 glue_1.3.1
## [16] reshape2_1.4.3 forcats_0.4.0 stringr_1.4.0
## [19] dplyr_0.8.3 purrr_0.3.3 readr_1.3.1
## [22] tidyr_1.0.0 tibble_2.1.3 ggplot2_3.2.1
## [25] tidyverse_1.3.0
##
## loaded via a namespace (and not attached):
## [1] readxl_1.3.1 backports_1.1.5
## [3] Hmisc_4.2-0 plyr_1.8.5
## [5] lazyeval_0.2.2 splines_3.6.0
## [7] crosstalk_1.0.0 BiocParallel_1.20.1
## [9] GenomeInfoDb_1.22.0 robust_0.4-18.2
## [11] digest_0.6.23 foreach_1.4.7
## [13] htmltools_0.4.0 GO.db_3.10.0
## [15] gdata_2.18.0 fansi_0.4.1
## [17] magrittr_1.5 checkmate_1.9.4
## [19] memoise_1.1.0 fit.models_0.5-14
## [21] cluster_2.0.8 doParallel_1.0.15
## [23] annotate_1.64.0 modelr_0.1.5
## [25] matrixStats_0.55.0 colorspace_1.4-1
## [27] blob_1.2.0 rvest_0.3.5
## [29] rrcov_1.4-7 haven_2.2.0
## [31] xfun_0.8 crayon_1.3.4
## [33] RCurl_1.95-4.12 jsonlite_1.6
## [35] genefilter_1.68.0 zeallot_0.1.0
## [37] impute_1.60.0 survival_2.44-1.1
## [39] iterators_1.0.12 gtable_0.3.0
## [41] zlibbioc_1.32.0 XVector_0.26.0
## [43] DelayedArray_0.12.2 BiocGenerics_0.32.0
## [45] DEoptimR_1.0-8 scales_1.1.0
## [47] mvtnorm_1.0-11 DBI_1.1.0
## [49] Rcpp_1.0.3 xtable_1.8-4
## [51] htmlTable_1.13.1 foreign_0.8-71
## [53] bit_1.1-15.1 preprocessCore_1.48.0
## [55] Formula_1.2-3 stats4_3.6.0
## [57] htmlwidgets_1.5.1 httr_1.4.1
## [59] ellipsis_0.3.0 acepack_1.4.1
## [61] farver_2.0.3 XML_3.98-1.20
## [63] pkgconfig_2.0.3 reshape_0.8.8
## [65] nnet_7.3-12 dbplyr_1.4.2
## [67] locfit_1.5-9.1 later_1.0.0
## [69] labeling_0.3 tidyselect_0.2.5
## [71] rlang_0.4.2 AnnotationDbi_1.48.0
## [73] munsell_0.5.0 cellranger_1.1.0
## [75] tools_3.6.0 cli_2.0.1
## [77] generics_0.0.2 RSQLite_2.2.0
## [79] broom_0.5.3 fastmap_1.0.1
## [81] evaluate_0.14 yaml_2.2.0
## [83] bit64_0.9-7 fs_1.3.1
## [85] robustbase_0.93-5 caTools_1.17.1.2
## [87] nlme_3.1-139 mime_0.8
## [89] xml2_1.2.2 compiler_3.6.0
## [91] rstudioapi_0.10 reprex_0.3.0
## [93] geneplotter_1.64.0 pcaPP_1.9-73
## [95] stringi_1.4.5 lattice_0.20-38
## [97] Matrix_1.2-17 vctrs_0.2.1
## [99] pillar_1.4.3 lifecycle_0.1.0
## [101] data.table_1.12.8 bitops_1.0-6
## [103] httpuv_1.5.2 GenomicRanges_1.38.0
## [105] R6_2.4.1 latticeExtra_0.6-28
## [107] promises_1.1.0 KernSmooth_2.23-15
## [109] IRanges_2.20.2 codetools_0.2-16
## [111] MASS_7.3-51.4 gtools_3.8.1
## [113] assertthat_0.2.1 SummarizedExperiment_1.16.1
## [115] DESeq2_1.26.0 withr_2.1.2
## [117] S4Vectors_0.24.2 GenomeInfoDbData_1.2.2
## [119] mgcv_1.8-28 parallel_3.6.0
## [121] hms_0.5.3 grid_3.6.0
## [123] rpart_4.1-15 rmarkdown_1.14
## [125] Cairo_1.5-10 shiny_1.4.0
## [127] Biobase_2.46.0 lubridate_1.7.4
## [129] base64enc_0.1-3